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I. INTRODUCTION 

Quantum computers are much more powerful than 
the classical ones 1,2 . However, a practical realization 
of such machines is still a big challenge ahead due to 
the fragility of the qubits and their decoherence arising 
from inevitable interactions with the environment. 
Over the years, several error correcting schemes have 
been developed 3-5 to make the proper working of a 
quantum computer possible, despite such interactions . 
Unfortunately, error correcting methods are themselves 
error prone due to the imperfections of the very gate 
which are responsible for detecting and correcting errors. 
Even this can be taken into account by the clever 
techniques of fault-tolerant quantum computation 5 . 
Unfortunately the error threshold below which fault 
tolerant quantum computation is possible, is very low'. 
In order to overcome such problems, one can combine 
the main quantum feature of the quantum world, namely 
superposition of states, with the robustness of classical 
bits which is the result of a macroscopic number of 
very small entities, comprising each bit. This has led 
to the idea of topological quantum computation 8,9 , 
where qubits are formed from degenerate ground states 
of a topologically ordered many body system. The 
important property that these degenerate ground states 
cannot be distinguished by local measurements, protects 
them agains local errors arising from coupling with 
the environment, at least as long as temperature is 
sufficiently low (see below). 

The prototype of such topological qubits is the toric 
code first introduced by Kitaev 8 . Since then other 
models have also been developed, notably the color code 
of Bombin and Martin-Delgado 10 . In these models, 
information is stored in the topologically degenerate 
ground states of the system and the computation is 
performed by braiding the quasiparticles (QP) of the 
model. Here one can use the gapped ground state and 
the robust nature of the topological phase of the model 


to protect information against local errors. The only 
perturbations that cause logical error are those with 
length equal to the system size. It should be noted that 
this applies only to zero temperatures, since it has been 
shown in 11 and elaborated in 12-14 that the proliferation 
of topological defects is capable of destroying topological 
order and hence quantum information encoded in such 
order at any finite temperature. Of course by tuning the 
parameters of such models 12 it is possible to increase 
the time scale for keeping the quantum information if 
the temperature is kept sufficiently low. 

In his seminal paper, Kitaev showed that universal 
quantum computation is not possible with an Abelian 
group and in order to perform all universal gates of 
quantum computation, one needs to either resort to the 
non-Abelian case or use other procedures like magic 
state distillation 15,16 or add some other non-topological 
resources like measurements to Abelian models as 
in 1 ' -20 . Aside from the ability of Abelian models to 
perform universal quantum computation in the above 
sense, their robustness and stability against external 
perturbations is still a crucial question which has to be 
investigated. Plenty of recent studies have been devoted 
for investigating such questions 21-31 . The motivation 
for such studies are not entirely based on quantum 
computation. In fact since topological order adds an 
entirely new paradigm for studying phase transition in 
condensed matter compared with the traditional one 
of symmetry breaking, it is very instructive to study 
various facets of this new phenomena. For example 
how topological order is destroyed under thermal 
fluctuation 11-13 , how it gives way to local order in the 
presence of an external field like magnetic field 2 ' are 
just two such broad questions. Of particular interest to 
us is the question of how topological order gives way to 
local order in the absence of external field when a local 
ordering interaction is added. This question was first 
investigated in our previous work 21 , where we studied 
the competition of topological order and ferromagnetic 
order in a model which we called Kitaev-Ising model. 
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There it was shown that there is a critical coupling of 
Ising where topological order is replaced with ferromag¬ 
netic order. In that model we used a mean-field analysis 
and a simple calculation of Wilson-Loop operators to 
find the transition point. 

In the context of condensed matter physics, it is 
important to understand the role of local degrees of 
freedom and their symmetry in topological order and 
its properties. It is mainly for such reasons, aside 
from the appeal of d-level states or qudits for quantum 
computation, that investigation of the above questions 
for general topological models for d — level states is im¬ 
portant. For example it is very desirable to understand 
the competition of Z d topological order and Z d local 
ferromagnetic order. The latter is a phase possessed 
by the d— level Potts model where the Zd symmetry 
is broken. The Potts model is a direct generalization 
of Ising model and has been extensively studied since 
its inception 3 '. While it is very desirable to perform 
such a study for general values of d , the analysis turns 
out to be very difficult. The 3—State Potts model is 
an exception in the sense that the Hamiltonian is much 
simpler, because it has an equidistant spectrum which 
brings much simplification in our analysis. It is therefore 
natural that Z 3 Potts model has been studied more 
intensively in the series of d— level Potts model, both in 
the classical context and in the new context in relation 
to topological order. 

It has been shown that the Z 3 Kitaev model is 
more robust against temperature than the Z 2 one 30 . 
The robustness of the model has also been studied in 
transverse magnetic field 31 and it has been shown that 
the topological order transforms to the magnetic order 
via a first order phase transition. Here we want to 
extend these studies by studying the competition of Z 3 
topological order induced by the Z 3 Kitaev model and 
the local order induced by the 3-State Potts model. 
We will see that the equidistant spectrum of the model 
leads to much technical simplification in the method we 
use, namely Perturbative Continuous Unitary Trans¬ 
formations (PCUT) 32-35 . We will show that like the 
Kitaev-Ising model 21 , the Kitaev-Potts model also shows 
a first-order phase transition and the topological order 
break downs at a critical coupling of Potts. In contrast 
to the Kitaev-Ising model, here our study is more 
comprehensive in that we use three different methods 
for locating the transition point, namely we capture 
the phase transition by analysis of the ground state 
energy, one-quasiparticle gap and geometric measure of 
entanglement. The location of the transition turns out 
to be slightly different by these methods, and we provide 
arguments to show that with better approximations of 
these methods, this slight change will be removed. 

The outline of the paper is as follows: In Sec. II and 
III, we briefly review the Z^ Kitaev and Potts models and 


their essential features needed for our study. In Sec. IV, 
we show that the low energy sector of the Z3 Kitaev- 
Potts model can be mapped to the 3-State Potts model 
in transverse magnetic field. We start our analysis of the 
mapped model in Sec. V by using mean held approxi¬ 
mation and Perturbative Continuous Unitary Transfor¬ 
mations method. By applying the PCUT method to the 
small Potts coupling limits, we calculate the one-quasi 
particle (1-QP) gap and also the Geometric Measure of 
Entanglement. We further compute the ground-state en¬ 
ergy in both small and large Potts couplings and capture 
the phase transition and the breakdown of the topologi¬ 
cal phase of the Kitaev model by analysis of the ground 
state energy and its derivatives, as well as the Geomet¬ 
ric Measure of Entanglement and energy gap. Finally, 
Sec. VI is devoted to the conclusion. 


II. Z d KITAEV MODEL 

7Ld Kitaev model is the generalization of the Kitaev 
model from Z 2 to Z d group 36 . The model can be defined 
on any two-dimensional lattice. In the present work, we 
consider the model on a square lattice on which the pe¬ 
riodic boundary conditions are imposed on both sides, 
where the lattice becomes a torus. The lattice has N 
plaquettes, N vertices and 2TV edges. The qudits live on 
the edges of the lattice. The Hamiltonian of the model 
consists of two kinds of operators, i.e. the plaquette 
and vertex operators. The operators are define based on 
the generalized Pauli operators (acting on a qudit) as: 
cr x \j) = \j + 1 mod d) and a z \j) = u 3 \j), w = e 27rz/d . 
These operators are not Hermitian and they don’t square 
to I except for d=2. They further obey the following 
commutation relation: a z a x = wo x tj z . The Hamiltonian 
of the model is given by the sum of the plaquette and 
vertex operators as: 

dd Kitaev • — j^. + 4)-^E( b p + b p), (!) 

S p 

where s and p denote the stars (vertices) and plaquettes 
respectively. In order for the model to be exactly solv¬ 
able, the A s ’s and B p ’s are defined such that they com¬ 
mute with each other. To this end an arbitrary direction 
is assigned to each edge and the A s and B p operators are 
defined as follows: 

• A s := Jlies a x\t ^th e link’s direction is inward, a x 
is applied, otherwise a~ x is, 

• b p : = ru P crz \ , by starting from a link and mov¬ 
ing counterclockwise, if each link’s direction is the 
same as moving’s direction cj z is applied , otherwise 
ctJ 1 is. (Fig. 1) 

These conventions lead to the commutativity of all A s ’s 
and B p s. Let us note that these arbitrary orientations do 
not have a physical significance, i.e. they lead to unitary 
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FIG. 1. (Color Online)The conventions for defining vertex 
and plaquette operators. A a := Yliea a xl’ where the + (—) 
corresponds to links with inward (outward) direction and 
B p := where the + (—) corresponds to counter¬ 

clockwise (clockwise) direction of the links around plaquettes. 


equivalent models. That is if we change the orientations 
on some links ( L ) in an arbitrary manner, we will end up 
with a model (H) which is iso-spectral with the original 
one ( H ). This unitary operator should transform a x to 
<j' x and <r z to in the Hamiltonian: 

H = UHU\ (7 = 0^, (2) 

ieL 

where 

Ui<T x ,iu\ = a f xi , UiO- z ,iuj = al ti , (3) 

It is straightforward to show that u is given by: 

U = Y1 \ d ~ k )( k \- ( 4 ) 

k 

Note that there are 2 N stabilizers in the Hamiltonian, 
but only 2 N — 2 of them are independent, because of the 
following two constraints on the torus: 

Y[A s = l[B p = I. (5) 

s p 

So there are d 2 degenerate ground states. The ground 
state is the state that is stabilized by all of the star and 
plaquette operators simultaneously and is equal to: 

1 00 ) := [](1 +A S + A 2 S +A 3 S + ... + Hf-^IO)® 2 ^. ( 6 ) 

S 

In order to construct the other d 2 — 1 degenerate ground 
states, we define the following four string operators: 

• 7),. ;| = n ieCi a z 1 > ^ starting from a link on the 
loop Ci and moving on it, if each link’s direction 
is the same as moving’s direction a z is applied , 
otherwise a~ x is, 

• t z , 2 = n ieC 2 a zi' by starting from a link on the 
loop C 2 and moving on it, if each link’s direction 
is the same as moving’s direction <r z is applied , 
otherwise a~ x is, 


C 2 



FIG. 2. (Color Online)Four nontrivial loops in a torus. Ci 
and O 2 are in the real lattice, while Ci and C 2 are in the dual 
lattice. 


• t xA = n ieCi a xl ’ ^ the links’s direction is the 
same as moving’s direction on Ci, a x is applied, 
otherwise a~ x is, 

• t x ,2 = n ieCo a xl ’ if the links’s direction is the 
same as moving’s direction on C 2 , o x is applied, 
otherwise a~ x is, 


where Ci, C 2 , Ci and C 2 are shown in Fig. 2. All the 
degenerate ground states are given as: 


1 ij) 


i,j = 0,1,2,..., d — 1. 

(7) 

One can further check that: 




II 

T z ,2\ij) = u J \ ij). 

(8) 


III. POTTS MODEL 


Considering any lattice of interest, the classical Potts 
model is defined by the following Hamiltonian 3 ': 

Hpotts = 'y \ d Si,sj j (9) 

id 


where the sum runs over the nearest neighbor sites of 
the lattice (i,j) and Si takes d different values. One can 
take these values to be d different roots of unity. For 
d= 2, Si = {1, —1}, which reduces the Hamiltonian (9) to 
the renowned Ising model. This simplification is a conse¬ 
quence of the following definition for the delta function: 
$Si,Sj = |(1 + SiSj). Such a definition can be extended 
to d-level and the analogue formula for the Potts model 
reads: 


d—l 


K,, = Td Z^r + (’Wn 


( 10 ) 


r—0 
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where s* is the complex conjugate of Sj. For general 
value of d , the model posses Zd symmetry s 3 —> uj.Sj. 

The quantum Potts model consists of d-level spins (qu- 
dits). The nearest neighbor spins interact with each other 
by the following Hamiltonian: 

d— 1 

Hpotts = — 7 ^ ^ ( a z,i a z,j) r ) ■ (11) 

(i,j) r=0 

The ground state of the system is the state where all 
spins are polarized in the same direction. So there are d 
degenerate ground states as follows: 

|<) := \i)® L , * = 0 , 1 , d—1, ( 12 ) 

where L is the number of spins. 

IV. THE 3-STATE KITAEV-POTTS 

Our aim is to study the phase transition of the Kitaev 
model in presence of the Potts interaction. As our se¬ 
ries expansion technique, i.e. the PCUT method is only 
applicable to those models with equidistant spectrum, 
we restrict our study to the d = 3 or qutrits and show 
that Z 3 Kitaev and 3-State Potts models have equidis¬ 
tant spectrum. The full Hamiltonian of the Z 3 Kitaev 
model perturbed by the Potts interaction is given by: 

H = HKitaev T ^Hpotts 

= -jJ2( A s + A\) - A^(B p + Bt) 

S p 

- ^ (( a zd a l,j) r + ( vm ^ t )- ( 1 3 ) 

(id) r=0 

where A is the perturbation parameter and is a mea¬ 
sure of the strength of the Potts interaction. The per¬ 
turbed Kitaev Hamiltonian is no longer exactly solvable. 
This is due to the fact that the a z operators in the 
Potts model do not commute with vertex operators of 
the Kitaev model. However, the plaquette operators still 
commute with the full Hamiltonian (13). The ground 
state of the Hamiltonian (13) is therefore in the sector in 
which B p = 1 for all plaquette operators. In this sector 
the Kitaev-Potts’s Hamiltonian reduces to the following 
form: 

H = -JY j { A s+A\) 

s 

- £ H (( a ^tdY + ( a z,i a zd) r ) - 2KN ■ ( 14 ) 

(id) r=0 

In order to tackle the Hamiltonian (14), we first define 
the following new basis and rewrite the full Hamiltonian 
in this new basis: 

|r) = |n,r 2 , -,r N ) ■= J|A’' ! '|0)® 2Ar , n = 0,1,2, 
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FIG. 3. (Color Online)In the mapped Kitaev-Potts Hamilto¬ 
nian, the original lattice (the lattices with -\/N * \fN Plaque- 
ttes, wherein y/~N is even) is de-coupled into two sublattices 
shown with blue dots and red dashes. 


where i runs over all the vertices of the lattice. This 
transformation was first defined in 21 for the the Kitaev- 
Ising model. One can readily check that B p \ r) = |r), 
We can therefore interpret Eq. (15) as the basis for the 
sector in which B p = 1 for all p. Next, we determine the 
action of the different terms in the Hamiltonian (14) on 
this new basis: 

A s |r) = A s J|A[*|O) 02Ar = ki,r 2 ,...,r g + l,...,rjv), 

(H3) 

so A s in this basis acts like generalized Pauli operator, cr x , 
which we denote by A'. Clearly, the action of Aj can be 
regarded as X'. In order to recast the Potts interaction 
in the new basis, we attach an arrow to each link. The 
direction of the arrows are illustrated in Fig. 2. The 
action of a Potts interaction term like a z ^ a o I b in the new 
basis, is denoted by: 

**,a< 6 |r )=<T Zia v\ ib A?A?Ag n ^riO> 02Ar (17) 

*^1,2,3 

= w ri w- r 3 |r)> 

which means that the Potts interaction term, a z ^ a <j\ b , 
commutes with all of the vertex operators except A\ and 
A 3 , actually it acts like ZiZ\ in this new basis: 

a z,a<^t ib = Z iZl. (18) 

It is straightforward to check that the following relations 
further hold in the new basis: 

<T t,a <T z,b = Z{Z 3 , 

vl, c az,d = Z\ Z 3 ■ 

a z,c a t,d ~ ^1^3- (19) 

Therefore the Z 3 Kitaev-Potts model in the new basis is 
given by: 

H = H a + H b - 2 AW, 


(15) 


( 20 ) 
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where 

h a = -j^(x 1+ x})-± y, E ( M t ) r + (4^ r , 

iSA (ij')£lr=0 

( 21 ) 

and 

^ = - J E(^+E)-^ E DM) r + (M) r - 

ieB (i,j)eB r=0 

( 22 ) 

The Hamiltonian (20) is nothing but the sum of two 
Potts models in a transverse magnetic field, wherein the 
Potts interactions act on nearest neighbor vertices in the 
two de-coupled sublattices shown in Fig. 3. 

The two de-coupled Hamiltonians are exactly the same. 
In the following sections we present our results for Ha- 
The results can be extended to the full Hamiltonian (20), 
without loss of generality. From now on to the end of the 
paper by H, we mean Ha, unless stated otherwise. 


V. METHODS 

In this section, we present the solution to the 
mapped model by applying the Mean field approxima¬ 
tion and Perturbative Continuous Unitary Transforma¬ 
tions method to the Hamiltonian (21). 


A. Mean Field approximation 

Suppose we are interested in finding the ground state 
energy of a given Hamiltonian. We therefore need to 
minimize (T|H|^) over all |\P)’s in the Hilbert space. 
One approximation (Mean field approximation) would be 
that we do the minimization only over product states. 
But the Hamiltonian has translational symmetry so we 
can search for the minimum energy only over states of the 
form In') = |$)® n , which have translational symmetry. 
For the mapped Hamiltonian n is the number of vertices 
in each sublattice, n= y. 

E = min® n ($|H|$)® n , (23) 

|$) 

where |$) is a general one-qutrit state. For the mapped 
Hamiltonian, 


H = -jY(Xi+Xj) 

i 


^ED^]) r +( 4 't. 

{i,j) r=0 


we have: 

E 0 = -Jn((X) + {X'))-4n±Y((Z r )(Z~ r ))- ( 24 ) 



FIG. 4. (Color Online) Ground state energy per site, eo, as 
a function of x. The inset represents its first derivative. The 
sharp jump in the derivation of energy at x c ~ 0.115 signals 
the first order quantum phase transition. 


Taking |<f>) = X)i=o a *N)> we should minimize: 



(25) 

It’s easy to check that the above equation recasts into: 


Eq — —Jn^|oo+ai+a 2 |‘" — 1^ —4,n\ ^|ao| 4 +|ai| 4 + |u 2 | 4 ^ ; 

(26) 

that has a permutation symmetry by changing the a^’s 
indices from (0,1,2) to (cr(l), a(2), <r(3)). If the state 
that minimize this expression respects the full symmetry, 
it should be equal to |<E>) = ^(|0) + |1) + |2)) and the 
ground state energy per site is equal to: 

£ 0 = -2J-^X. (27) 

Setting J = 1, we see that Eq. (27) is analytic for any 
value of A and one does not observe any phase transition. 
However, if some of the symmetries are broken and the 
ground state remains invariant under only one permuta¬ 
tion, it will be in the form of |$) = a 0 |0) + ai|l) + ai|2). 
By setting ao = sin0, ai = ^e* a cos0, the ground state 
energy per site reduces to: 

e 0 = —J{ cos 2 9 + \pl sin 29 cos a) — 4A(sin 4 9 + ^ cos 4 9). 

(28) 

By minimizing eg, the ground state as a function of per¬ 
turbation parameter x = (this explicit form of x is 
chosen because it makes the comparison between mean- 
field and PCUT results easier) is given as follows: 


|cf>) 


^(|!) + | 2 » 

sin0*|O) + ^|l> + ^|2> 


for x < 0.115 


for x > 0.115 
(29) 
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where, 

sin 29 X - 2y/2 cos 2 9 X 

x = --—^-=-- —(30) 

36 sin 29 x (sin 2 9 X — \ cos 2 9 X ) 

The above relation shows that the nature of the ground 
state changes at a critical point x c = 0.115, which can 
be a signal of the phase transition. Finally we analyze 
the case where the symmetry is fully broken and the 
ground state has the general form as |$) = XEo 0 *!*)- 
In this case, we minimize the energy numerically. Figure 
4 demonstrates the mean-field ground state energy per 
site as a function of x = |j. As we can see, there is 
a small kink in the ground state energy curve (a sharp 
jump in its derivative) at x c ~ 0.115. So the results 
show that the ground state actually respects some per¬ 
mutation symmetries. This jump signals a first order 
quantum phase transition. In order to investigate the 
critical point more accurately, we resort to a more ac¬ 
curate approximation technique, i.e. the Perturbative 
Continuous Unitary Transformation method 32-35 which 
is the subject of the next subsection. 

B. Perturbative Continuous Unitary 
Transformations 

In this section, we briefly review the Perturbative Con¬ 
tinuous Unitary Transformation (PCUT) method and 
apply it to the small- and large-coupling limits of the 
problem. 

The Continuous Unitary Transformation (CUT) 
method which was first introduced by Wegner 32 in the 
framework of condensed matter theory, is basically used 
to diagonalize or block-diagonalize a given Hamiltonian 
by applying an infinite number of unitary operators to 
the initial Hamiltonian in a continuous fashion as: 

H{t) = U\£)HU{£), (31) 

where t is the continuous flow parameter such that H = 
H(£ = 0) and H 0 g = H(£ = oo) is the (block-) diago¬ 
nal Hamiltonian. The Hamiltonian is transformed by a 
unitary operator which its evolution is governed by: 

diU{i) = (32) 

In which rj{t) is the anti-Hermitian generator of the uni¬ 
tary transformation U{t). Combining Eq. (31) and (32) 
together one can show that the initial Hamiltonian flows 
in the form of a differential commutator equation : 

d e H(e) = [ V (£),H(l)\. (33) 

The method therefore requires the choice of a suit¬ 
able generator for the unitary operators to obtain the 
desirable form of H e g. Ulirig and Knetter introduced 
the quasiparticle (QP) conserving generator which is very 


well suited for our purpose. We refer the interested reader 
to Refs 33,34 for detailed discussions on QP conserving 
generator. 

The perturbative version of the CUT method (PCUT), 
can be applied to the Hamiltonians of the form H = 
Q + XV where the first part of the Hamiltonian, Q 1 is 
diagonal with an equidistant spectrum bounded from be¬ 
low and the second part can be treated as a perturbation 
(A is the expansion parameter). The method further re¬ 
quires that the perturbing part can be written in the 
form V = En=-N T n , where T n increments (decrements, 
if n < 0) the number of excitations (quasiparticles) by 
n such that [Q,T n ] = nT n 33 . Transforming the initial 
problem by using the QP conserving generator, the effec¬ 
tive Hamiltonian is brought to the form that conserves 
the number of quasiparticles, [H e s,Q] = 0. The energy 
spectrum of the system can therefore be extracted per- 
turbatively by acting the H e g on the ground state and 
the multi-particle sectors of the Hilbert space. 

1. Small-coupling limit (X <C J) 

In the following we discuss the procedure of applying 
the PCUT method to the small-coupling limit (A -C J) 
of the mapped Kitaev-Potts model, i.e. the Potts model 
in transverse magnetic field on the A or B sublattices of 
Fig. 3: 

H = -jY / (X i + Xj) 

i£A 

E E((^jr + (Mr)- ( 34 ) 

(i,j)eA r =0 

The first term in Eq. (34) is an effective field term which 
is diagonal. Denoting the local vacuum of each site by 
10), the elementary excitations of the model for d = 3 
are two separate spin flips labeled by | 1 ) and | 2 ) which 
correspond to the eigenstates of X operator with w and 
w -1 eigenvalues, respectively. Let us note that either of 
the excitations cost an energy of 3 J. The elementary ex¬ 
citations are energetically indistinguishable. So the first 
term for d = 3 has an equidistant spectrum and can 
be regarded as Q for implementation of the PCUT, and 
the PCUT results obtained in 1-QP sector of the Hilbert 
space such as 1-QP gap are degenerate for both of the 
excitations. 

The Potts interaction at the right side of Eq. (34) can 
be treated as a perturbation V which for d = 3 is denoted 
by: 

F = -y E (&%) + (M)) + (35) 

(i,j) 

where C is a constant. The perturbing part consists of 
two-body interactions and can change the number of ex¬ 
citations over the ground state of the effective field term 
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Q by n = {0, ±1, ±2} when it acts on the bonds of the 
square lattice. Therefore the Hamiltonian (34) can be 
written as: 

H = Q ~ x(T 2 +Ti + T 0 + T_i + T_ 2 ), (36) 

where 

<3 = £ -<*■++) (37) 

i 

is the quasiparticle counting operator (/ is the identity 
operator), x = is the expansion parameter and T n 
operators are given by: 



r +2 = |12) iJ <00| +121)^(001, (38) 

(id) 

T+i = E | 22 ) iJ < 01 | + 122)^(101 

(id) 

+ 111)^(021 + 111 )^( 201 , 

To = 53 112)^(211 + | 01 ) M -( 10 | + | 02 )jj( 20 | 

(id) 

+ 121)^(121 + 110)^(011 + 120 )^( 021 . 

From the hermiticity condition = T_„, and absorb¬ 
ing the expansion parameter — x in the definition of T n 
operators, Eq. (36) is recast into: 

n—2 

H = Q + 53 Tn - ( 3g ) 

n=— 2 

Under the CUTs, the above Hamiltonian is continuously 
transformed with the flow parameter t as: 

n—2 

H(£) = Q+ 53 T n {t), (40) 

n=— 2 

where we wish to reach a situation where T n (t = oo) = 0 
for all n / 0. In order to fulfil this demand, we choose 
the quasiparticle conserving generator in the following 
form 33 ’ 34 : 

VW = T+ 2 (£) + T +1 (£) - T_ x (£) - T_ 2 (£). (41) 

With this choice of generator, the flow Eq. (33) can be 
written as: 

d e T 0 (£) = 2 [T +2 {£),T. 2 {£)} + 2[T + 1 (£),T_ 1 (£)], (42) 

d t T +2 {l) = -2 T +2 {£) + {T +2 (£),T 0 (£)\, 

d t T +1 {£) = -T +1 (£) + 2 [T +2 (£),T. 1 (£)\ + [T +1 (£),T 0 (£]}. 

Let us stress that the Hamiltonian remains hermition 
under the unitary transformation. We can therefore cal¬ 
culate T_ 2 (£) and T_i(£) from the hermiticity condition. 
Solving the flow equation is still a very cumbersome task 
because there are an infinite number of terms in T n {£). 


FIG. 5. (Color Online) 1-QP dispersion of the Kitaev-Potts 
model in the small coupling regime. The minimum of the 
dispersion is located at K=(0,0). 

However, we can tackle the problem by performing a per¬ 
turbative expansion of the flow equation. We can there¬ 
fore introduce the expansion of the T n operators as 35 : 

OO 

T„(£)=53TW(£), (43) 

i= 1 

where i is the order of perturbation. Using this relation, 
the perturbative expansion of the flow equation is written 
as: 

dt T 0 (fe) (£) = 2 53 [T$ (£), T%- j) (£)] + 2 ^ [T^ (£), (£)], 

j=i j =i 

deT™(£) = -2 T<$(£) + ^[T^(£),T 0 (fc - j) (£)], 

3 = 1 
fc-1 

d t T™ (£) = -T$ {£) + 2 53 [T+2 (£), (£)] ( 44 ) 

j=i 

+ j2iTl j he)X k - j) m- 

3 = 1 

Solving this set of equations for the initial condition 
Tn\i = 0) = 8\^T n and then taking the limit of (£ —» 
oo), we can obtain H c g. The effective QP conserving 
Hamiltonian for the small-coupling limit up to order 3 in 
perturbation parameter is obtained as: 

2 

= Q — xT 0 + x 2 ^, T_ x ] + y [T 2 , T_ 2 ] 

-y (W I t 0’T-2]} + [[r 2 ,r 0 ],r_ 2 ]) 

-y (pi, P~i, T_ 2 ]] + (45) 

— ([TulTo^ + iiTuTolT^y 
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FIG. 6. (Color Online) The 1-QP gap of the Kitaev-Potts 
model for different DlogPade approximants. The critical point 
highlighted by gray dashed line is consistent with the closure 
of the gap at x c ~ 0.129. The inset demonstrates the scaling 
on the closure of the bare series of the gap for different orders 
of perturbation. 


The ground state energy and 1-QP gap can be obtained 
by acting the effective hamiltonian on OP and IP sector 
of the Q. We have calculated the ground state energy 
and 1-QP energy gap of the system in the small-coupling 
limit up to order 8 in perturbation theory: 
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(46) 


A sc = 1 - 4a; - 10a; 2 - 5a; 3 - 


1895a; 4 14107a; 5 


6 


3572759a; 6 26566267a; 7 


216 


324 


+ 18 
85919559673a 8 
77760 


(47) 


Fig. 6 illustrates the 1-QP gap of the system in the 
small-coupling limit as a function of a = |j. The bare 
and extrapolated series are well converged. Closure of 
the 1-QP gap occurs at a c ~ 0.129. The inset of Fig. 6 
further demonstrates scaling of the closure of the bare 
series of the gap for different orders of perturbation. As 
one can see, by increasing the order of perturbation the 
transformation point starts to become smaller until it 
converges to x c ~ 0.12 which is fully consistent with the 
closure of DlogPad extrapolants. 

Existence of anyonic excitations in the system is one 
of the signatures of the topological order. As one can see 
form Fig. 6, the anyone gap is open until in the vicinity of 
the transition point. It is therefore reasonable to point 
out that the topological order survives during the tun¬ 
ing of the perturbation until the critical point at which 
the anyons condenses due to the strength of the Potts 
interaction and the topological order breaks down 38,39 . 

Let us further note that the knowledge of the gap is 
not solely sufficient to determine the first- or second- 
order nature of the phase transition 38,39 and one has to 


FIG. 7. (Color Online) The ground state energy per site £o 
as a function of 9. Solid line corresponds to the bare series of 
order 8. Symbols denote different Pade approximants. 


analyze the ground state energy of the system and its 
derivatives 40 to capture the phase transition correctly. 
We therefore postpone further discussion on the phase 
transition to Sec. VB3 after we calculate the ground 
state energy series in large-coupling limits. 


2. Large-coupling limit (\^> J) 


Following our discussion, we now apply the PCUT 
method to the large-coupling limit of the problem (A 
J). For J = 0, the Hamiltonian is 3-State Potts model 
which is ferromagnetically ordered and has an equidis¬ 
tant spectrum thus satisfying the first condition for the 
PCUT method. The excitations of the model further 
correspond to the anti-ferromagnetic bonds of the square 
lattice. When J / 0, the effective field term in Eq. (34) 
can be considered as a perturbation which changes the 
number of excitations by n = {0, ±1, ±2, ±3, ±4}. We 
can therefore write the Hamiltonian of the model in the 
large-coupling limit in terms of T n operators as: 

4 

H = Q~hJ2 T n, (48) 

n =—4 

where Q is the quasiparticle counting operator defined 
as: 


M-'Eii4)(Z i zh + {Z j Zl) 

Q = - ( *' 3> J — * , (49) 

o 

and h = J/2X is the expansion parameter. Using the 
PCUT method, we have calculated the ground state en¬ 
ergy per site, £q, in the large-coupling limit up to order 
8 in perturbation parameter by acting the H e g- on the 
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FIG. 8. (Color Online) First derivative of the ground state 
energy per site ( d\Eo) as a function of 6. Solid line corre¬ 
sponds to the bare series of order 8. Symbols denote different 
Pade approximants. 
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FIG. 9. (Color Online) Second derivative of the ground state 
energy per site ( d\eo) as a function of 9. Solid line corre¬ 
sponds to the bare series of order 8. Symbols denote different 
Pade approximants. 


0-QP sector of the Hilbert space: 


8 h 2 h 3 19h 4 3/i 5 12779/i 6 

3 2 8 672 128 846720~ ' 

1052987 h 7 458808396457/1® 

121927680 ~~ 62768369664000' 


As we have already mentioned in the previous section, 
one can obtain a better understanding about the nature 
of the phase transition by analysis of the ground state 
energy results. We will therefore provide strong evidences 
for the first-order phase transition in the Z 3 Kitaev-Potts 
model in the next section. 


3. Analyzing series expansion results 

To investigate the nature of the phase transition, we 
study the ground state energy per site of the system in 
the whole range of the expansion parameter by merging 
the small- and large-coupling results. Setting A = sin 8 
and J = cos 9 , we can join the small- and large-coupling 
results to obtain a complete picture. Fig. 7 demonstrates 
the ground state energy per site £0 as a function of 8. 
The small- and large-coupling series cross each other at 
9 C ss 0.52, giving rise to a kink in the £0 curve which is 
fully consistent with the first-order phase transition. The 
location of the kink is essentially the same for different 
Pade extrapolants and has strong agreement with the 
closure of the 1-QP gap ( x c = | tan$ c ss 0.129). 

Using the Feynman-Hellman theorem, we have cal¬ 
culated the first derivative of £ 0 which is equivalent to 
the magnetization in statistical mechanics. Figs. 8 and 
9 depict first and second derivatives of £0 for different 
9 values. The sharp jumps in the first derivative is 
a clear signature of the first-order phase transition 40 . 
Let us further note that the first order nature of the 


phase transition in the Z 3 Potts model in transverse mag¬ 
netic field has already been confirmed using series ex¬ 
pansion combined with infinite projected entangled-pair 
state (IPEPS) method in Ref. 31 . The jump in the second 
derivative further confirms the location of the transition 
point which is exactly the same as that of the ground- 
state energy for the bare series and different Pade extrap¬ 
olants. 


C. Geometric Measure of Entanglement 

In this section we calculate a measure of multipar¬ 
tite entanglement, Geometric Measure of Entanglement 
(GME) 41 , as another tool for capturing the phase tran¬ 
sition. GME is a measure of multipartite entanglement 
in quantum many body systems which measures the dis¬ 
tance, Hilbert Schmidt distance, between a given state, 
^(x)), and the closest product state, |P), as follows: 

GME = -log 2 ^max|(P|'F(x))| 2 ^ , (51) 

where the maximization is over all product states. This 
implies that the more entangled the states are, more dis¬ 
tance from the set of product states, which is a convex 
set, they will have. In order to calculate the GME, we 
have to calculate the ground state of the mapped Hamil¬ 
tonian, |\P(cc)), by PCUT procedure up to a specific or¬ 
der of perturbation. The PCUT method transforms the 
initial Hamiltonian by a unitary transformation ( [U{oo )) 
to a block-diagonal form (H c g), the basis of the Hilbert 
space is rotated, such that H 0 g commutes with Q. So 
they have the same eigenstates. It’s also known that 
H e g’s ground state is the vacuum state of Q with no ex¬ 
citations (|0)) 42 . The H 0 g and the initial Hamiltonian, 
H , are unitary equivalent and their ground states are 
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related to each other as follows: 


jvhOr)} = [/(oc)|0). (52) 

In order to determine U(f), we expand Eq. (32) as a 
function of the perturbation parameter x: 

k-1 

diU (k \i) = -^[/ w (^ (h) (4 (53) 

j =o 

where 

v (k-j) (q = T cfc-i) (£) + yg-i) (£) _ rg _i) (£) - 2 -J) (£), 

(54) 

and k is the order of perturbation. By using Eqs. (44), 
and the initial condition U k { 0) = lJ^o, one can solve the 
above equation perturbatively up to second order in x. 
Taking the limit i —> oo, the U(o o) is found to be: 


u = 1 + |(T +2 - t_ 2 ) + x(t +1 - r_o 

+ ^(T +1 -T_ 1 ) 2 + y(T +2 -T_ 2 ) 2 

+ ^[To, (T +2 + r_ 2 )] + x 2 [T 0 , (T +1 + T_i)] 

+ ^(T +1 T +2 + T— \ T— 2 + 2 T +2 T +1 + 2T_ 2 T_ 1 ) 

6 

+ y([Ti,T_ 2 ] + [T_!,T 2 ] - T_ 2 Ti - T 2 T_!). (55) 


The ground state of the initial Hamiltonian H is there¬ 
fore given by: 


|f(x)) = t/|0) = (l- 4 ^)|6) + ^|l 2) 


X 

2 


4 n 


X 

4 


El 1 2 ) + 2 E 


4n 


8 n 


1 o 
0 2 


4n 


T 


4n 


2 2 
2 0 


2 E 

4 n 


l l 

1 o 


+ 2^|2 2 2) + 2^|l 1 1) 


(56) 


■2EI 1 0 2 ) 


t( E 

\ 8n(2n—7) 


1 2 
1 2 


This state is a superposition of some product states, these 
product states are eigenstates of Q. The first term, |0) 
refers to the ground state of Q with no excitations, which 
is the state of all spins in the X eigenstate with +1 
eigenvalue, |0). |l 2) is the state of all the spins in 
|0) state, except for the two of them, which are nearest 
neighbour, and have the states of |1) and |2). The sum¬ 
mation must be done over all the states of this kind, and 
one should note that in a lattice with n vertices the num¬ 
ber of this kind of states are 2 times the number of bonds, 



FIG. 10. (Color Online)GME as a function of x for n=25. 
Convexity of GME changes close to the x c ~ 0.16. The inset 
further demonstrates the derivative of GME as a function of 
x. The sharp jump in the figure is a clear signature of phase 
transition. 


i.e. 2 x 2n = 4n. Furthermore, 


1 0 
0 2 


and |l 0 2) 


denote the state with all spins in the |0) state, except for 
the two of them, which are next nearest neighbours and 


have the states of |1) and |2). The last term, 


1 2 
1 2 


also refers to the state with all spins in |0) state except 
for the four of them. The state of these four spins consist 
of two clusters of 11 2), where these two clusters can be 

separate or not. 


IT (a;)) is exact up to the second order in perturba¬ 
tion parameter x and is normalized. We have also calcu¬ 
lated the ground state of the initial Hamiltonian up to 
the fourth order, numerically and calculated the GME. 
Since the ground state has translational symmetry, the 
closest product state also preserves this symmetry 43 . 
We can therefore perform the maximization only over 
the states in the form | P) = |^)® n , where \(j>) = 
cos(0)|O)+e _l “ sin(0) sin(c/j)|l) + e - *^ sin(0) cos(<p)|2). In 
other word, the maximization is only over 4 parameters 
which makes the numerical calculation of the GME pos¬ 
sible. Fig. 10 illustrates the geometric measure of entan¬ 
glement as a function oi x. As we can see, the convexity 
of the GME changes sign from positive to negative close 
to a critical point, x c ~ 0.16. The sharp jump in the 
derivative of GME captures the phase transition more 
clearly, since it’s known that for two-dimensional systems 
discontinuity in the derivative of multipartite entangle¬ 
ment leads to a quantum phase transition 44 . The critical 
point is very close to the one obtained from analysis of 
the ground state energy and the gap i.e. x c ~ 0.129. The 
slight difference is likely the consequence of the difference 
in orders of perturbation. This can be best inferred by 
reminding the scaling behavior of the closure of 1-QP gap 
(inset of Fig. 6) which predicts that increasing the order 
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of perturbation makes the result to become smaller and 
converge to ss 0.12. The same fact should also holds here 
for the location of sharp jump in the derivative of the 
GME. 


VI. CONCLUSION 

The Z d Kitaev model is a system with a topologically 
ordered ground state which, when complemented with 
some other resources like magic state distillation 15,16 
or measurements 20 , is well suited for the purpose of 
universal quantum computation without resorting to 
non-Abelian groups. It is therefore of great interest to 
study the stability and robustness of the topological 
phase of the model in the presence of external per¬ 
turbations. In this paper, we studied the competition 
between the topological order induced by the Z 3 Kitaev 
model and the local order induced by the 3-State Potts 
model on the square lattice. We showed that the 
Kitaev model in the presence of the Potts interaction 


is mapped to Potts model in a transverse magnetic 
field. Using the high-order series expansion based on the 
continuous unitary transformations in the small- and 
large- Potts couplings, we showed that the topological 
phase breaks down to a non-topological phase with 
local order parameter through a first-order quantum 
phase transition at x c = ~ 0.129. Our results 

were further in good agreement with the mean-field 
approximation results at x c ss 0.115. Computing 
the Geometric Measure of Entanglement also shows 
that the derivative of GME has a sharp jump very 
close to the critical point and also, the convexity of 
GME changes sign from positive to negative in this point. 
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